ProgettoLampedusa

# Carica le librerie necessarie
library(here)         # Gestione dei percorsi dei file
library(corrplot)     # Visualizzazione matrice di correlazione
library(nortest)      # Test di normalità
library(lmtest)       # Test diagnostici per modelli lineari
library(car)          # Test di collinearità e altri test per modelli lineari
library(plotly)       # Visualizzazioni interattive
library(heatmaply)    # Heatmap interattivo
library(ggheatmap)    # Heatmap con ggplot2
library(gridExtra)    # Organizzazione di grafici con grid
library(ggfortify)    # Visualizzazioni grafiche per modelli statistici
library(olsrr)        # Diagnostica e grafici per modelli di regressione
library(ggplot2)      # Creazione di grafici con ggplot2
library(viridis)      # Colormaps viridis
library(glmnet)       # Modelli di regressione con penalizzazione lasso ed elastic net
library(highcharter)  # Creazione di grafici interattivi con Highcharts
library(reshape2)     # Manipolazione dati
library(RColorBrewer) # Colormaps
library(tidyr)        # Manipolazione dati
library(dplyr)        # Manipolazione dati

# Configura il percorso di root per la generazione di report
knitr::opts_knit$set(root.dir = here("0_Materiale"))

INIZIALIZZAZIONE DATI E GRAFICI DATI

# Leggi il dataset da un file di testo
dataset <- read.delim("basketball_teams.txt")

# Definisci il primo e l'ultimo anno del range da considerare per lo studio
FIRST <- 1976
LAST <- 2011

# Filtra il dataset secondo le condizioni specificate
df <- dataset[dataset$lgID == "NBA" & dataset$year >= FIRST & dataset$year <= LAST & dataset$games == 82,]

# Converti la colonna lgID in un fattore per consentire la generazione di variabili dummy
dataset$lgID <- as.factor(dataset$lgID)

# Stampiamo un riassunto statistico del dataframe filtrato
summary(df)

##       year          lgID               tmID             franchID        
##  Min.   :1976   Length:841         Length:841         Length:841        
##  1st Qu.:1985   Class :character   Class :character   Class :character  
##  Median :1993   Mode  :character   Mode  :character   Mode  :character  
##  Mean   :1993                                                           
##  3rd Qu.:2001                                                           
##  Max.   :2008                                                           
##     confID             divID                rank          confRank     
##  Length:841         Length:841         Min.   :0.000   Min.   : 1.000  
##  Class :character   Class :character   1st Qu.:2.000   1st Qu.: 4.000  
##  Mode  :character   Mode  :character   Median :3.000   Median : 7.000  
##                                        Mean   :3.565   Mean   : 7.164  
##                                        3rd Qu.:5.000   3rd Qu.:10.000  
##                                        Max.   :8.000   Max.   :15.000  
##    playoff              name               o_fgm          o_fga     
##  Length:841         Length:841         Min.   :2565   Min.   :5972  
##  Class :character   Class :character   1st Qu.:2981   1st Qu.:6592  
##  Mode  :character   Mode  :character   Median :3220   Median :6903  
##                                        Mean   :3239   Mean   :6941  
##                                        3rd Qu.:3489   3rd Qu.:7253  
##                                        Max.   :3980   Max.   :8868  
##      o_ftm          o_fta          o_3pm           o_3pa            o_oreb    
##  Min.   :1189   Min.   :1475   Min.   :  0.0   Min.   :   0.0   Min.   : 720  
##  1st Qu.:1523   1st Qu.:2039   1st Qu.: 78.0   1st Qu.: 293.0   1st Qu.: 974  
##  Median :1649   Median :2201   Median :283.0   Median : 814.0   Median :1083  
##  Mean   :1666   Mean   :2212   Mean   :279.1   Mean   : 804.6   Mean   :1086  
##  3rd Qu.:1797   3rd Qu.:2371   3rd Qu.:445.0   3rd Qu.:1267.0   3rd Qu.:1191  
##  Max.   :2388   Max.   :3051   Max.   :837.0   Max.   :2283.0   Max.   :1520  
##      o_dreb         o_reb          o_asts          o_pf          o_stl       
##  Min.   :2044   Min.   :2922   Min.   :1422   Min.   :1476   Min.   : 455.0  
##  1st Qu.:2348   1st Qu.:3381   1st Qu.:1760   1st Qu.:1777   1st Qu.: 613.0  
##  Median :2433   Median :3506   Median :1934   Median :1900   Median : 671.0  
##  Mean   :2434   Mean   :3520   Mean   :1934   Mean   :1909   Mean   : 681.3  
##  3rd Qu.:2527   3rd Qu.:3647   3rd Qu.:2094   3rd Qu.:2033   3rd Qu.: 746.0  
##  Max.   :2966   Max.   :4216   Max.   :2575   Max.   :2470   Max.   :1059.0  
##       o_to          o_blk           o_pts           d_fgm          d_fga     
##  Min.   : 910   Min.   :204.0   Min.   : 6901   Min.   :2488   Min.   :5638  
##  1st Qu.:1207   1st Qu.:360.0   1st Qu.: 7958   1st Qu.:2978   1st Qu.:6593  
##  Median :1311   Median :411.0   Median : 8404   Median :3243   Median :6911  
##  Mean   :1338   Mean   :421.8   Mean   : 8423   Mean   :3239   Mean   :6941  
##  3rd Qu.:1443   3rd Qu.:473.0   3rd Qu.: 8879   3rd Qu.:3493   3rd Qu.:7268  
##  Max.   :2011   Max.   :716.0   Max.   :10371   Max.   :4265   Max.   :8142  
##      d_ftm          d_fta           d_3pm            d_3pa          d_oreb    
##  Min.   :1217   Min.   :  0.0   Min.   :   0.0   Min.   :1579   Min.   : 745  
##  1st Qu.:1514   1st Qu.: 82.0   1st Qu.: 282.0   1st Qu.:2033   1st Qu.: 986  
##  Median :1648   Median :280.0   Median : 836.0   Median :2203   Median :1095  
##  Mean   :1666   Mean   :279.1   Mean   : 804.6   Mean   :2212   Mean   :1086  
##  3rd Qu.:1808   3rd Qu.:446.0   3rd Qu.:1251.0   3rd Qu.:2395   3rd Qu.:1177  
##  Max.   :2377   Max.   :683.0   Max.   :1768.0   Max.   :3071   Max.   :1495  
##      d_dreb         d_reb          d_asts          d_pf          d_stl      
##  Min.   :2012   Min.   :2976   Min.   :1336   Min.   :1434   Min.   :461.0  
##  1st Qu.:2326   1st Qu.:3378   1st Qu.:1778   1st Qu.:1788   1st Qu.:623.0  
##  Median :2431   Median :3497   Median :1939   Median :1900   Median :677.0  
##  Mean   :2434   Mean   :3520   Mean   :1934   Mean   :1909   Mean   :681.3  
##  3rd Qu.:2529   3rd Qu.:3653   3rd Qu.:2092   3rd Qu.:2020   3rd Qu.:734.0  
##  Max.   :3067   Max.   :4309   Max.   :2537   Max.   :2453   Max.   :955.0  
##       d_to          d_blk           d_pts        o_tmRebound  d_tmRebound
##  Min.   : 949   Min.   :264.0   Min.   : 6909   Min.   :0    Min.   :0   
##  1st Qu.:1208   1st Qu.:380.0   1st Qu.: 7968   1st Qu.:0    1st Qu.:0   
##  Median :1304   Median :419.0   Median : 8453   Median :0    Median :0   
##  Mean   :1338   Mean   :421.8   Mean   : 8423   Mean   :0    Mean   :0   
##  3rd Qu.:1444   3rd Qu.:460.0   3rd Qu.: 8841   3rd Qu.:0    3rd Qu.:0   
##  Max.   :1980   Max.   :654.0   Max.   :10723   Max.   :0    Max.   :0   
##     homeWon         homeLost        awayWon         awayLost        neutWon 
##  Min.   : 6.00   Min.   : 1.00   Min.   : 1.00   Min.   : 8.00   Min.   :0  
##  1st Qu.:21.00   1st Qu.:10.00   1st Qu.:10.00   1st Qu.:21.00   1st Qu.:0  
##  Median :26.00   Median :15.00   Median :15.00   Median :26.00   Median :0  
##  Mean   :25.63   Mean   :15.37   Mean   :15.37   Mean   :25.63   Mean   :0  
##  3rd Qu.:31.00   3rd Qu.:20.00   3rd Qu.:20.00   3rd Qu.:31.00   3rd Qu.:0  
##  Max.   :40.00   Max.   :35.00   Max.   :33.00   Max.   :40.00   Max.   :0  
##     neutLoss    confWon         confLoss         divWon         divLoss     
##  Min.   :0   Min.   : 5.00   Min.   : 7.00   Min.   : 1.00   Min.   : 1.00  
##  1st Qu.:0   1st Qu.:20.00   1st Qu.:20.00   1st Qu.: 8.00   1st Qu.: 9.00  
##  Median :0   Median :27.00   Median :26.00   Median :12.00   Median :12.00  
##  Mean   :0   Mean   :26.91   Mean   :26.91   Mean   :12.27   Mean   :12.27  
##  3rd Qu.:0   3rd Qu.:34.00   3rd Qu.:33.00   3rd Qu.:16.00   3rd Qu.:16.00  
##  Max.   :0   Max.   :48.00   Max.   :52.00   Max.   :25.00   Max.   :27.00  
##       pace             won          lost        games         min       
##  Min.   :  0.00   Min.   :11   Min.   :10   Min.   :82   Min.   :19680  
##  1st Qu.:  0.00   1st Qu.:31   1st Qu.:32   1st Qu.:82   1st Qu.:19780  
##  Median :  0.00   Median :42   Median :40   Median :82   Median :19805  
##  Mean   :  6.71   Mean   :41   Mean   :41   Mean   :82   Mean   :19817  
##  3rd Qu.:  0.00   3rd Qu.:50   3rd Qu.:51   3rd Qu.:82   3rd Qu.:19855  
##  Max.   :102.00   Max.   :72   Max.   :71   Max.   :82   Max.   :20080  
##     arena             attendance       bbtmID         
##  Length:841         Min.   :    0   Length:841        
##  Class :character   1st Qu.:32767   Class :character  
##  Mode  :character   Median :32767   Mode  :character  
##                     Mean   :32728                     
##                     3rd Qu.:32767                     
##                     Max.   :32767

# STUDIO DATI RELATIVI AI RIMBALZI

# o_oreb = Rimbalzi ottenuti in attacco
# o_dreb = Rimbalzi subiti in attacco
# o_reb  = totale rimbalzi in attacco
# d_oreb = Rimbalzi subiti in difesa
# d_dreb = Rimbalzi ottenuti in difesa
# d_reb  = totale rimbalzi in difesa

# Seleziona le colonne pertinenti dal dataframe
df_reb <- subset(df, select = c("o_oreb", "o_dreb", "o_reb", "d_oreb", "d_dreb", "d_reb", "won"))

# Imposta il layout a 2 righe e 3 colonne per i grafici di densità
par(mfrow = c(2, 3))

# Ciclo per generare grafici di densità per ciascuna variabile di interesse
for (variables in 1:(dim(df_reb)[2]-1)) {
  thisvar <- df_reb[, variables]
  d <- density(thisvar)
  xmin <- floor(min(thisvar))
  xmax <- ceiling(max(thisvar))
  
  # Crea il plot della densità con stile accattivante
  plot(d, main = names(df_reb)[variables], xlab = "", col = "blue", lwd = 1.5, xlim = c(xmin, xmax), ylim = c(0, max(d$y)*1.1))
  
  # Aggiungi la distribuzione normale teorica ideale in rosso
  x <- seq(xmin, xmax, length = 100)
  lines(x, dnorm(x, mean = mean(thisvar), sd = sd(thisvar)), col = "red", lwd = 1.5)  # Modifica lo spessore delle linee
  
  # Aggiungi griglia per migliorare la leggibilità
  grid()
}

# Aggiungi titolo in grassetto e corsivo, con spessore del testo modificato
title(bquote(bold(italic("Density plots with Normal Distribution"))), line = -17, cex.main = 2, outer = TRUE)

# Stampa percentuale di valori non zero per ciascuna variabile di interesse
print(paste("Percentage non-zero o_oreb: ", round(length(which(df_reb$o_oreb > 0)) / dim(df_reb)[1] * 100, 2)))

## [1] "Percentage non-zero o_oreb:  100"

print(paste("Percentage non-zero o_dreb: ", round(length(which(df_reb$o_dreb > 0)) / dim(df_reb)[1] * 100, 2)))

## [1] "Percentage non-zero o_dreb:  100"

print(paste("Percentage non-zero o_reb: ", round(length(which(df_reb$o_reb > 0)) / dim(df_reb)[1] * 100, 2)))

## [1] "Percentage non-zero o_reb:  100"

print(paste("Percentage non-zero d_oreb: ", round(length(which(df_reb$d_oreb > 0)) / dim(df_reb)[1] * 100, 2)))

## [1] "Percentage non-zero d_oreb:  100"

print(paste("Percentage non-zero d_dreb: ", round(length(which(df_reb$d_dreb > 0)) / dim(df_reb)[1] * 100, 2)))

## [1] "Percentage non-zero d_dreb:  100"

print(paste("Percentage non-zero d_reb: ", round(length(which(df_reb$d_reb > 0)) / dim(df_reb)[1] * 100, 2)))

## [1] "Percentage non-zero d_reb:  100"

STUDIO SULLE VITTORIE

SANDBOX TESTING

### SANDBOX

# Creazione di un vettore per archiviare un elenco di squadre
teams <- c()

# Creazione di una nuova variabile 'reb' come somma di 'o_reb' e 'd_reb'
df$reb <- c(df$o_reb + df$d_reb)

# Visualizzazione di un riepilogo statistico sulla variabile 'reb'
summary(df$reb)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    5963    6805    7002    7040    7245    8359

# Loop per ottenere un elenco unico di squadre
for (team in df$tmID) {
  teams <- unique(c(teams, team))
}

# Ordinamento delle squadre
teams <- sort(teams)

# Calcolo di valori aggregati per le variabili specificate
values <- aggregate(cbind(reb, o_oreb, o_dreb, d_oreb, d_dreb, o_reb, d_reb, won) ~ tmID, data = df, FUN = sum)

# Creazione di dati per una distribuzione normale
x <- seq(-4, 4, by = 0.01)
y <- dnorm(x, mean = 0, sd = 1)
normal <- data.frame(x, y)

# ISTOGRAMMA INTERATTIVO
histogram <- plot_ly(data = df, x = ~reb, type = "histogram",
        marker = list(color = 'skyblue', line = list(color = 'black', width = 1)),
        nbinsx = 20, legendgroup = "Rimbalzi", name = "Campione") %>% 
  add_trace(type = "scatter", mode = "lines",
            x = c(mean(df$reb), mean(df$reb)),
            y = c(0, 215),
            line = list(color = "red", width = 2, dash = "dash"),
            name = "Media") %>%
  layout(title = list(text = "<b><i>Istogramma dei Rimbalzi</i></b>", pad = 10),
         xaxis = list(title = '<b><i>Rimbalzi</i></b>'),
         yaxis = list(title = '<b><i>Frequenza</i></b>'),
         legend = list(title = "<b><i>Legenda</i></b>", tracegrouporder = "reversed"))

(histogram)

## A marker object has been specified, but markers is not in the mode
## Adding markers to the mode...

## Warning: 'scatter' objects don't have these attributes: 'nbinsx'
## Valid attributes include:
## 'cliponaxis', 'connectgaps', 'customdata', 'customdatasrc', 'dx', 'dy', 'error_x', 'error_y', 'fill', 'fillcolor', 'fillpattern', 'groupnorm', 'hoverinfo', 'hoverinfosrc', 'hoverlabel', 'hoveron', 'hovertemplate', 'hovertemplatesrc', 'hovertext', 'hovertextsrc', 'ids', 'idssrc', 'legendgroup', 'legendgrouptitle', 'legendrank', 'line', 'marker', 'meta', 'metasrc', 'mode', 'name', 'opacity', 'orientation', 'selected', 'selectedpoints', 'showlegend', 'stackgaps', 'stackgroup', 'stream', 'text', 'textfont', 'textposition', 'textpositionsrc', 'textsrc', 'texttemplate', 'texttemplatesrc', 'transforms', 'type', 'uid', 'uirevision', 'unselected', 'visible', 'x', 'x0', 'xaxis', 'xcalendar', 'xhoverformat', 'xperiod', 'xperiod0', 'xperiodalignment', 'xsrc', 'y', 'y0', 'yaxis', 'ycalendar', 'yhoverformat', 'yperiod', 'yperiod0', 'yperiodalignment', 'ysrc', 'key', 'set', 'frame', 'transforms', '_isNestedKey', '_isSimpleKey', '_isGraticule', '_bbox'

# PLOT DI DENSITA CON SOVRAPPOSIZIONE DI UNA NORMALE IDEALE
density_data <- density(df$reb)
mu <- mean(df$reb)
sigma <- sd(df$reb)
normal_data <- dnorm(density_data$x, mean = mu, sd = sigma)

density_plot <- plot_ly(x = density_data$x, y = density_data$y, type = 'scatter', mode = 'lines',
             line = list(color = 'blue', width = 2),
             name = "Densità rimbalzi totale") %>%
  layout(title = "Densità dei Rimbalzi",
         xaxis = list(title = "Rimbalzi"),
         yaxis = list(title = "Density", autotick = TRUE, autorange = TRUE))

density_plot <- add_trace(density_plot, x = density_data$x, y = normal_data, mode = 'lines',
               line = list(color = 'green', width = 2),
               fill = "tozeroy", fillcolor = "rgba(0, 255, 0, 0.2)",
               name = "Dist normale ideale")

density_plot <- add_trace(density_plot, x = c(mu, mu), y = c(0, max(density_data$y)),
               mode = 'lines', line = list(color = 'red', width = 2, dash = 'dash'),
               name = "Media")

(density_plot)

# CORRPLOT
data_subset <- df[, c("reb", "o_reb", "d_reb")]
cor_matrix <- cor(data_subset)

rownames(cor_matrix) <- colnames(data_subset)
colnames(cor_matrix) <- colnames(data_subset)

cor_data <- reshape2::melt(cor_matrix)
names(cor_data) <- c("Var1", "Var2", "Corr")

dimnames(cor_matrix) <- list(rownames(cor_matrix), colnames(cor_matrix))
data_for_plotly <- as.data.frame(as.table(cor_matrix))

cor_plot <- plot_ly(data = cor_data, 
             x = ~Var1, 
             y = ~Var2, 
             z = ~Corr, 
             type = "heatmap", 
             colors = colorRampPalette(c("#4575b4", "#91bfdb", "#e0f3f8", "#fee08b", "#d73027"))(100),
             hoverinfo = "x+y+z") %>% 
  layout(title = 'Correlation Matrix',
         xaxis = list(title = "", tickangle = 45, side = "bottom", automargin = TRUE),
         yaxis = list(title = "", automargin = TRUE),
         autosize = TRUE)

cor_values <- round(as.matrix(cor_matrix), 2)  # Round for readability
for (i in seq_len(nrow(cor_matrix))) {
  for (j in seq_len(ncol(cor_matrix))) {
    cor_plot <- cor_plot %>% add_annotations(
      x = rownames(cor_matrix)[i],
      y = colnames(cor_matrix)[j],
      text = as.character(cor_values[i, j]),
      showarrow = FALSE,
      font = list(color = ifelse(cor_values[i, j] < 0.5, "white", "black"))
    )
  }
}
(cor_plot)

# BOX PLOT
aggregated_data <- df %>%
  select(tmID, year, o_reb, d_reb) %>%
  group_by(tmID, year) %>%
  summarize(
    o_reb = mean(o_reb, na.rm = TRUE), 
    d_reb = mean(d_reb, na.rm = TRUE),
    .groups = "drop"  # Aggiunto per evitare il raggruppamento
  )
reshaped_data <- aggregated_data %>%
  pivot_longer(cols = c(o_reb, d_reb), names_to = "stat_type", values_to = "stat") %>%
  unite("new_col", tmID, stat_type, sep = "_") %>%
  pivot_wider(names_from = new_col, values_from = "stat")
box_plot <- plot_ly()
for (i in 2:ncol(reshaped_data)) {
    current_column_data <- reshaped_data[[i]]
    box_plot <- box_plot %>% add_trace(y = current_column_data, name = colnames(reshaped_data)[i], type = "box")
}

(box_plot)

## Warning: Ignoring 30 observations

## Warning: Ignoring 30 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 19 observations

## Warning: Ignoring 19 observations

## Warning: Ignoring 4 observations

## Warning: Ignoring 4 observations

## Warning: Ignoring 23 observations

## Warning: Ignoring 23 observations

## Warning: Ignoring 8 observations

## Warning: Ignoring 8 observations

## Warning: Ignoring 24 observations

## Warning: Ignoring 24 observations

## Warning: Ignoring 12 observations

## Warning: Ignoring 12 observations

## Warning: Ignoring 13 observations

## Warning: Ignoring 13 observations

## Warning: Ignoring 1 observations

## Warning: Ignoring 1 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 29 observations

## Warning: Ignoring 29 observations

## Warning: Ignoring 30 observations

## Warning: Ignoring 30 observations

## Warning: Ignoring 31 observations

## Warning: Ignoring 31 observations

## Warning: Ignoring 31 observations

## Warning: Ignoring 31 observations

## Warning: Ignoring 13 observations

## Warning: Ignoring 13 observations

## Warning: Ignoring 9 observations

## Warning: Ignoring 9 observations

## Warning: Ignoring 26 observations

## Warning: Ignoring 26 observations

## Warning: Ignoring 1 observations

## Warning: Ignoring 1 observations

## Warning: Ignoring 19 observations

## Warning: Ignoring 19 observations

## Warning: Ignoring 3 observations

## Warning: Ignoring 3 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 21 observations

## Warning: Ignoring 21 observations

## Warning: Ignoring 11 observations

## Warning: Ignoring 11 observations

# HEATMAP
heatmap_df <- subset(values, select = -c(o_oreb, o_dreb, d_oreb, d_dreb, won, reb))
rownames(heatmap_df) <- heatmap_df$tmID
heatmap_df <- heatmap_df[,-1]

heatmap_plot <- heatmaply(
 heatmap_df,
 colors = viridis(n = 256,  option = "magma"),
 k_col = 2,
 k_row = 4,
)

(heatmap_plot)

# BARPLOT
bar_plot <- plot_ly(values, x = ~tmID, y = ~o_reb, type = 'bar', name = 'Rimbalzi offensivi', marker = list(color = '#FFAFA1')) %>%
  add_trace(y = ~d_reb, name = 'Rimbalzi difensivi', marker = list(color = '#b2fff8')) %>%
  layout(yaxis = list(title = 'Valori'), barmode = 'stack')

(bar_plot)

#TEST ANDERSON-DARLING

# Questo codice esegue il test di Anderson-Darling sulla variabile 'reb' nel tuo dataset (df)

ad.test(df$reb)

## 
##  Anderson-Darling normality test
## 
## data:  df$reb
## A = 3.6997, p-value = 3.1e-09

#Con un livello di significatività (α) di 0.01 e un p-value molto piccolo (3.1e-09) ottenuto dal test di normalità di Anderson-Darling per i dati della variabile df\(reb, puoi concludere che hai sufficiente evidenza statistica per respingere lipotesi nulla che i dati seguono una distribuzione normale.Con il tuo livello di significatività del 0.01 e il p-value molto piccolo (3.1e-09), il p-value è inferiore al livello di significatività, quindi respingeresti lipotesi nulla. Questo suggerisce che i dati nella variabile df\)reb non seguono una distribuzione normale al livello di significatività del 0.01. In termini più pratici, hai abbastanza evidenza statistica per concludere che la variabile df$reb non segue una distribuzione normale basandoti sui risultati del test di Anderson-Darling.

#TEST KOLMOGOROV SMIRNOV

# Test di Kolmogorov-Smirnov per confrontare la distribuzione di 'reb' con una distribuzione normale
ks.test(df$reb, "pnorm")

## Warning in ks.test.default(df$reb, "pnorm"): i legami non dovrebbero essere
## presenti per il test di Kolmogorov-Smirnov

## 
##  Asymptotic one-sample Kolmogorov-Smirnov test
## 
## data:  df$reb
## D = 1, p-value < 2.2e-16
## alternative hypothesis: two-sided

Il risultato che hai ottenuto riguarda il test di Kolmogorov-Smirnov a campione singolo sui dati contenuti nella variabile df$reb. Il test KS confronta la distribuzione empirica dei tuoi dati con una distribuzione teorica (spesso una distribuzione uniforme). In breve, il risultato suggerisce che i tuoi dati non seguono la distribuzione teorica presunta, e cè un elevata probabilità che la differenza osservata sia statisticamente significativa.

#TEST SHAPIRO WILK

# Test di Shapiro-Wilk per la normalità dei dati nella variabile 'reb' (rimbalzi)
sf.test(df$reb)

## 
##  Shapiro-Francia normality test
## 
## data:  df$reb
## W = 0.98016, p-value = 1.758e-08

#In sintesi, il risultato del test di Shapiro-Francia indica che i tuoi dati nella variabile df$reb non seguono una distribuzione normale. Questo è supportato dal valore basso del p-value, il quale suggerisce che la differenza tra la distribuzione dei tuoi dati e una distribuzione normale è statisticamente significativa.

#STUDIO SULLE VITTORIE


# ISTOGRAMMA INTERATTIVO
# Calcola la media della variabile "won" nel data frame "df"
mean_value <- mean(df$won)
# Crea un istogramma con plot_ly
histogram <- plot_ly(df, x = ~won, type = "histogram",
                     marker = list(color = "skyblue", line = list(color = "white", width = 0.5)), 
                     opacity = 0.7, name = "Campione") %>%
  layout(title = list(text = "<b><i>Distribuzione delle Vittorie</i></b>", y = 0.97),
         xaxis = list(title = "<b><i>Numero di Vittorie</i></b>", zeroline = FALSE),
         yaxis = list(title = "<b><i>Frequenza</i></b>", zeroline = FALSE),
         barmode = "overlay") %>%
  add_trace(x = ~mean(won), type = "scatter", mode = "lines", 
            line = list(color = "red", width = 2), name = "<i><b>Media</i></b>") %>%
  add_trace(x = ~mean(won), type = "scatter", mode = "markers", 
            marker = list(color = "red", size = 8), showlegend = FALSE) %>%
  add_annotations(text = sprintf("<i><b>Media: %.2f</b></i>", mean(df$won)), x = mean_value, y = 0, 
                  arrowhead = 2, arrowcolor = "red", arrowsize = 1.5, arrowwidth = 2)

# Visualizza l'istogramma interattivo
histogram

## A marker object has been specified, but markers is not in the mode
## Adding markers to the mode...

# PLOT DENSITA INTERATTIVO

# Calcola la densità delle vittorie
density_plot <- density(df$won)

# Crea il plot di densità con plot_ly
density_interactive <- plot_ly(x = density_plot$x, y = density_plot$y, type = "scatter", mode = "lines",
                               line = list(color = "skyblue", width = 2), name = "Densità") %>%
  layout(title = list(text = "<b><i>Distribuzione di Densità delle Vittorie</i></b>", x = 0.5),
         xaxis = list(title = "<i>Numero di Vittorie</i>"),
         yaxis = list(title = "<i>Densità</i>"),
         showlegend = TRUE) %>%
  add_trace(x = c(mean(df$won), mean(df$won)), y = c(0, max(density_plot$y)),
            type = "scatter", mode = "lines", line = list(color = "red", width = 2),
            name = "<i><b>Media</b></i>") %>%
  add_trace(x = mean(df$won), y = max(density_plot$y), type = "scatter", mode = "markers",
            marker = list(color = "red", size = 8), showlegend = FALSE) %>%
  add_annotations(text = sprintf("<i><b>Media: %.2f</b></i>", mean(df$won)), x = mean(df$won), y = max(density_plot$y) * 1.05,
                  arrowhead = 2, arrowcolor = "red", arrowsize = 1.5, arrowwidth = 2,
                  ax = 0, ay = -40)

# Visualizza il plot di densità interattivo
density_interactive

## A line object has been specified, but lines is not in the mode
## Adding lines to the mode...

# CORRPLOT

data_subset <- df[, c(11:40, 54)]

cor_matrix <- cor(data_subset, use = "complete.obs")

cor_data <- melt(cor_matrix)
names(cor_data) <- c("Variable1", "Variable2", "Correlation")

# Crea la heatmap
p <- plot_ly(data = cor_data, x = ~Variable1, y = ~Variable2, z = ~Correlation,
             type = "heatmap",
             colors = viridis(n = 1024, option = "magma"),
             hoverinfo = "x+y+z") %>%
  layout(
    title = '<b><i>Correlation Matrix Heatmap</i></b>',
    xaxis = list(
      title = list(text = "<b><i>Variabile 1</i></b>", standoff = 0),
      tickangle = 45,
      zeroline = FALSE
    ),
    yaxis = list(
      title = list(text = "<b><i>Variabile 2</i></b>", standoff = 0),
      tickangle = 45,
      zeroline = FALSE
    ),
    autosize = TRUE
  )

# Visualizza la heatmap
(p)

# BOXPLOT

aggregated_data <- df %>%
  select(tmID, year, won) %>%
  group_by(tmID, year) %>%
  summarize(
    won = mean(won, na.rm = TRUE),
  )

## `summarise()` has grouped output by 'tmID'. You can override using the
## `.groups` argument.

reshaped_data <- aggregated_data %>%
  pivot_longer(cols = c(won), names_to = "stat_type", values_to = "stat") %>%
  unite("new_col", tmID, stat_type, sep = "_") %>%
  pivot_wider(names_from = new_col, values_from = "stat")

fig <- plot_ly()

for (i in 2:ncol(reshaped_data)) {
  current_column_data <- reshaped_data[[i]]
  fig <- fig %>% add_trace(y = current_column_data, name = colnames(reshaped_data)[i], type = "box")
}

# Visualizza il boxplot interattivo
(fig)

## Warning: Ignoring 30 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 19 observations

## Warning: Ignoring 4 observations

## Warning: Ignoring 23 observations

## Warning: Ignoring 8 observations

## Warning: Ignoring 24 observations

## Warning: Ignoring 12 observations

## Warning: Ignoring 13 observations

## Warning: Ignoring 1 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 29 observations

## Warning: Ignoring 30 observations

## Warning: Ignoring 31 observations

## Warning: Ignoring 31 observations

## Warning: Ignoring 13 observations

## Warning: Ignoring 9 observations

## Warning: Ignoring 26 observations

## Warning: Ignoring 1 observations

## Warning: Ignoring 19 observations

## Warning: Ignoring 3 observations

## Warning: Ignoring 27 observations

## Warning: Ignoring 21 observations

## Warning: Ignoring 11 observations

df$reb <- df$o_reb + df$d_reb

INIZIALIZZAZIONE MODELLO DI REGRESSIONE LINEARE

L’IMPORTANZA DEI RIMBALZI

# Questo codice effettua alcune operazioni di preprocessing sui dati, come la creazione di nuove variabili (f1 a f10) e la divisione del dataset in set di training e test. Successivamente,      viene creato un modello lineare (linMod) e ne viene eseguito il resampling. Infine, viene creato un modello lineare con le covariate normalizzate (linModNormalized) e vengono generati         grafici di diagnostica per entrambi i modelli.


# o_oreb = Rimbalzi ottenuti in attacco
# o_dreb = Rimbalzi subiti in attacco
# o_reb  = totale rimbalzi in attacco = o_oreb + o_dreb
# d_oreb = Rimbalzi subiti in difesa
# d_dreb = Rimbalzi ottenuti in difesa
# d_reb  = totale rimbalzi in difesa

# Definizione di nuove variabili
df$f1 <- (df$o_oreb) / (df$o_fga - df$o_fgm)
df$f2 <- (df$d_dreb) / (df$d_fga - df$d_fgm)
df$f3 <- (df$o_oreb + 1.5 * df$d_dreb) / (df$o_dreb + 2 * df$d_oreb)
df$f4 <- (df$o_oreb - df$o_dreb) + 1.5 * (df$d_dreb - df$d_oreb)
df$f5 <- (df$d_oreb / df$d_to) / (df$o_dreb / df$o_to)
df$f6 <- (df$o_oreb + df$d_dreb) - (df$d_oreb - df$o_dreb)^2
df$f7 <- (df$f1 + df$f2)^2
df$f8 <- (df$o_oreb) / (df$o_reb)
df$f9 <- ((df$o_oreb) / (df$o_dreb))^2
df$f10 <- ((df$d_dreb) / (df$d_oreb))^2


# Divisione in Test e Train per evitare che il modello fitti troppo bene sui nostri dati
sample <- sample(c(TRUE, FALSE), nrow(df), replace=TRUE, prob=c(0.7, 0.3))

df = subset(df, select = c("f1", "f2", "f3", "f4", "f5", "f6", "f7", "f8", "f9", "f10", "won", "divID", "confID"))

train  <- df[sample, ]
test   <- df[!sample, ]

#ATTENZIONE: facendo il train sul valore delle variabili, questo significa che sono esse ad essere i nostri dati, non i rimbalzi in se.
#            importante non confondersi durante l'esposizione.


# Creazione del modello lineare
linMod <- lm(won ~ f1 + f2 + f3 + f4 + f5 + f6 + f7 + f8 + f9 + f10, data = train)
summary(linMod)

## 
## Call:
## lm(formula = won ~ f1 + f2 + f3 + f4 + f5 + f6 + f7 + f8 + f9 + 
##     f10, data = train)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.7696  -3.1466   0.0134   3.1314  13.3878 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.021e+03  8.602e+01  11.874  < 2e-16 ***
## f1          -6.650e+02  1.481e+02  -4.491 8.52e-06 ***
## f2          -1.055e+03  1.565e+02  -6.744 3.65e-11 ***
## f3          -3.961e+02  3.630e+01 -10.911  < 2e-16 ***
## f4           7.341e-02  7.724e-03   9.504  < 2e-16 ***
## f5          -1.819e+02  4.754e+00 -38.253  < 2e-16 ***
## f6           1.236e-05  1.484e-06   8.328 5.64e-16 ***
## f7           4.622e+02  7.823e+01   5.908 5.81e-09 ***
## f8          -2.537e+02  7.397e+01  -3.429 0.000647 ***
## f9           3.281e+01  2.236e+01   1.467 0.142912    
## f10          4.732e+00  9.798e-01   4.829 1.74e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.683 on 597 degrees of freedom
## Multiple R-squared:  0.8602, Adjusted R-squared:  0.8578 
## F-statistic: 367.2 on 10 and 597 DF,  p-value: < 2.2e-16

# Normalizziamo le covariate
train$f1_z <- scale(train$f1)
train$f2_z <- scale(train$f2)
train$f3_z <- scale(train$f3)
train$f4_z <- scale(train$f4)
train$f5_z <- scale(train$f5)
train$f6_z <- scale(train$f6)
train$f7_z <- scale(train$f7)
train$f8_z <- scale(train$f8)
train$f9_z <- scale(train$f9)
train$f10_z <- scale(train$f10)

# Normalizzo i dati di test
test$f1_z <- scale(test$f1)
test$f2_z <- scale(test$f2)
test$f3_z <- scale(test$f3)
test$f4_z <- scale(test$f4)
test$f5_z <- scale(test$f5)
test$f6_z <- scale(test$f6)
test$f7_z <- scale(test$f7)
test$f8_z <- scale(test$f8)
test$f9_z <- scale(test$f9)
test$f10_z <- scale(test$f10)

linModNormalized <- lm(won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + f7_z + f8_z + f9_z + f10_z, data = train)

# Imposta il layout della pagina
par(mfrow = c(2, 2))
# Grafico dei residui standardizzati
plot(linModNormalized, which = 1, col = "skyblue", pch = 16, main = "Residui Standardizzati")
abline(h = 0, col = "red", lty = 2)  # Aggiungi una linea orizzontale a zero nel grafico dei residui standardizzati
# Grafico dei livelli
plot(linModNormalized, which = 2, col = "lightgreen", pch = 16, main = "Grafico dei Livelli")
# Grafico della distribuzione dei residui
plot(linModNormalized, which = 3, col = "pink", pch = 16, main = "Distribuzione dei Residui")
# Grafico di Q-Q plot dei residui
plot(linModNormalized, which = 4, col = "orange", pch = 16, main = "Q-Q Plot dei Residui")

# Aggiungi altre personalizzazioni a piacere, come titoli, etichette degli assi, ecc.
# Ripristina il layout predefinito
par(mfrow = c(1, 1))

TEST SUL MODELLO DI REGRESSIONE LINEARE

TEST BREUSCH-PAGAN (Test di omoschedasticità)

# TEST SUL MODELLO DI REGRESSIONE LINEARE

# 1 Summary
summary(linModNormalized)

## 
## Call:
## lm(formula = won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + 
##     f7_z + f8_z + f9_z + f10_z, data = train)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.7696  -3.1466   0.0134   3.1314  13.3878 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  40.4030     0.1899 212.755  < 2e-16 ***
## f1_z        -21.0213     4.6809  -4.491 8.52e-06 ***
## f2_z        -50.2841     7.4564  -6.744 3.65e-11 ***
## f3_z        -25.4430     2.3318 -10.911  < 2e-16 ***
## f4_z         19.5615     2.0582   9.504  < 2e-16 ***
## f5_z        -12.8531     0.3360 -38.253  < 2e-16 ***
## f6_z          6.4973     0.7802   8.328 5.64e-16 ***
## f7_z         38.0074     6.4329   5.908 5.81e-09 ***
## f8_z         -8.2180     2.3964  -3.429 0.000647 ***
## f9_z          2.0413     1.3915   1.467 0.142912    
## f10_z         6.8481     1.4180   4.829 1.74e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.683 on 597 degrees of freedom
## Multiple R-squared:  0.8602, Adjusted R-squared:  0.8578 
## F-statistic: 367.2 on 10 and 597 DF,  p-value: < 2.2e-16

# 2 R-quadrato e R-quadrato Adattato
summary_linModNormalized <- summary(linModNormalized)
r_squared <- summary_linModNormalized$r.squared
cat("R-squared:", r_squared, "\n")

## R-squared: 0.8601693

n <- length(df$o_oreb)
k <- length(linModNormalized$coefficients) - 1
adjusted_r_squared <- 1 - ((1 - r_squared) * (n - 1) / (n - k - 1))
cat("Adjusted R-squared:", adjusted_r_squared, "\n")

## Adjusted R-squared: 0.9872881

# 2 Test Shapiro per valutare la normalità dei residui
shapiro.test(residuals(linModNormalized))

## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(linModNormalized)
## W = 0.99432, p-value = 0.02275

# 3 Test di omoschedasticità (Breusch-Pagan test) --> il risultato suggerisce omoschedasiticità
bptest(linModNormalized)

## 
##  studentized Breusch-Pagan test
## 
## data:  linModNormalized
## BP = 14.173, df = 10, p-value = 0.1652

# 4 Test di multicollinearità
car::vif(linModNormalized)

##        f1_z        f2_z        f3_z        f4_z        f5_z        f6_z 
##  606.573045 1539.138319  150.527628  117.266989    3.125347   16.850589 
##        f7_z        f8_z        f9_z       f10_z 
## 1145.588930  158.982557   53.601265   55.664014

#TROVO GLI OUTLAYERS

# Il codice implementa l'individuazione degli outliers considerando i residui che superano una soglia moltiplicata per la deviazione standard.

# Calcola i residui dal modello di regressione lineare normalizzato
residui <- residuals(linModNormalized)

# Definisci una soglia per gli outlier
soglia_outlier <- 2
outliers <- which(abs(residui) > soglia_outlier*sd(residui))

# Identifica gli outlier basati sulla soglia definita
outliers <- which(abs(residui) > soglia_outlier * sd(residui))
outliers

##  440  470  519  652  769  799  822  837  891  900  914  929  936  940  945  957 
##    1   23   62  162  239  260  277  290  328  334  344  357  361  365  369  379 
## 1026 1107 1128 1139 1205 1214 1224 1301 
##  405  462  480  489  533  538  545  603

#MODELLO SENZA OUTLAYERS

#Il codice rimuove gli outliers, ricrea il modello lineare normalizzato e produce grafici di diagnostica per entrambi i modelli.


# Rimuovi gli outliers e ricrea il modello lineare normalizzato
if (length(outliers) != 0) {
  train_1 <- train[-outliers,]
} else {
  train_1 <- train
}

# ATTENZIONE: facendo il train sul valore delle variabili, questo significa che sono esse ad essere i nostri dati, non i rimbalzi in sé. Importante non confondersi durante l'esposizione.


linModNormalized_1 <- lm(won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + f7_z + f8_z + f9_z + f10_z, data = train_1)

# Stampa i summary dei modelli
summary(linModNormalized)

## 
## Call:
## lm(formula = won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + 
##     f7_z + f8_z + f9_z + f10_z, data = train)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.7696  -3.1466   0.0134   3.1314  13.3878 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  40.4030     0.1899 212.755  < 2e-16 ***
## f1_z        -21.0213     4.6809  -4.491 8.52e-06 ***
## f2_z        -50.2841     7.4564  -6.744 3.65e-11 ***
## f3_z        -25.4430     2.3318 -10.911  < 2e-16 ***
## f4_z         19.5615     2.0582   9.504  < 2e-16 ***
## f5_z        -12.8531     0.3360 -38.253  < 2e-16 ***
## f6_z          6.4973     0.7802   8.328 5.64e-16 ***
## f7_z         38.0074     6.4329   5.908 5.81e-09 ***
## f8_z         -8.2180     2.3964  -3.429 0.000647 ***
## f9_z          2.0413     1.3915   1.467 0.142912    
## f10_z         6.8481     1.4180   4.829 1.74e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.683 on 597 degrees of freedom
## Multiple R-squared:  0.8602, Adjusted R-squared:  0.8578 
## F-statistic: 367.2 on 10 and 597 DF,  p-value: < 2.2e-16

summary(linModNormalized_1)

## 
## Call:
## lm(formula = won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + 
##     f7_z + f8_z + f9_z + f10_z, data = train_1)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -9.0156 -2.9560  0.1958  3.0363  9.3808 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  40.1581     0.1742 230.573  < 2e-16 ***
## f1_z        -18.9076     4.7109  -4.014 6.77e-05 ***
## f2_z        -48.9506     7.4952  -6.531 1.44e-10 ***
## f3_z        -27.7666     2.1909 -12.674  < 2e-16 ***
## f4_z         22.5459     1.9693  11.449  < 2e-16 ***
## f5_z        -13.0032     0.3099 -41.956  < 2e-16 ***
## f6_z          7.1655     0.7371   9.721  < 2e-16 ***
## f7_z         36.0699     6.4990   5.550 4.37e-08 ***
## f8_z         -9.8409     2.1973  -4.479 9.07e-06 ***
## f9_z          2.5699     1.2670   2.028    0.043 *  
## f10_z         7.4600     1.3101   5.694 1.98e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.206 on 573 degrees of freedom
## Multiple R-squared:  0.884,  Adjusted R-squared:  0.882 
## F-statistic: 436.7 on 10 and 573 DF,  p-value: < 2.2e-16

# Grafici
ols_plot_resid_lev(linModNormalized)

ols_plot_resid_lev(linModNormalized_1)

# Influence plot per il primo modello
influencePlot(linModNormalized, id = 5)

## Warning in applyDefaults(id, defaults = list(method = "noteworthy", n = 2, :
## unnamed id arguments, will be ignored

##       StudRes         Hat       CookD
## 769 -2.337140 0.140832493 0.080791911
## 822  2.567135 0.141797087 0.098069596
## 940  2.887305 0.007421224 0.005597557
## 957  3.037847 0.108873451 0.101105957

# Personalizzazione del primo grafico
par(mar=c(5, 5, 2, 2))
title(main = "Influence Plot - Modello 1", col.main = "blue", font.main = 4)

# Apri una nuova finestra grafica
dev.new()

# Influence plot per il secondo modello
influencePlot(linModNormalized_1, id = 5)

## Warning in applyDefaults(id, defaults = list(method = "noteworthy", n = 2, :
## unnamed id arguments, will be ignored

##          StudRes        Hat        CookD
## 454  -1.92755178 0.07523360 2.734929e-02
## 592   0.41907733 0.12016398 2.183703e-03
## 643   2.26269646 0.04316141 2.084517e-02
## 903  -0.02277664 0.13549723 7.404726e-06
## 907   2.25139453 0.01156714 5.354468e-03
## 1206  1.56841793 0.10113810 2.509849e-02

# Personalizzazione del secondo grafico
par(mar=c(5, 5, 2, 2))
# Aggiungi un titolo sopra al grafico
title(main = "Influence Plot - Modello 2", col.main = "blue", font.main = 4)

# METODO DI REGRESSIONE LASSO

# Il codice implementa il metodo di regressione LASSO con cross-validation per trovare il miglior valore lambda.

# Definizione della variabile di risposta
y <- train_1$won

# Definizione della matrice delle variabili predittive (uso solo poche variabili ma potete farlo con tutte da togliere però la risposta area)
x <- data.matrix(train_1[, c("f1_z", "f2_z", "f3_z", "f4_z", "f5_z", "f6_z", "f7_z", "f8_z", "f9_z", "f10_z")])

# Esegui la cross-validation k-fold per trovare il valore lambda ottimale
cv_model <- cv.glmnet(x, y, alpha = 1)

# Trova il valore lambda ottimale che minimizza il test MSE
best_lambda <- cv_model$lambda.min
best_lambda

## [1] 0.003616601

# Addestramento del modello con cross-validation
cv_model <- cv.glmnet(x, y)

# Grafico personalizzato
plot(cv_model, xvar="lambda", main="", xlab="log(Lambda)", ylab="Mean Squared Error", col="blue", lwd=2)

## Warning in plot.window(...): parametro grafico "xvar" non valido

## Warning in plot.xy(xy, type, ...): parametro grafico "xvar" non valido

## Warning in axis(side = side, at = at, labels = labels, ...): parametro grafico
## "xvar" non valido

## Warning in axis(side = side, at = at, labels = labels, ...): parametro grafico
## "xvar" non valido

## Warning in box(...): parametro grafico "xvar" non valido

## Warning in title(...): parametro grafico "xvar" non valido

# Aggiungi un titolo personalizzato con formattazione LaTeX per grassetto e corsivo
title(main=expression(bold(italic("Validazione Incrociata"))), line = 3)

# Fittiamo il modello con il miglior lambda (penalizzazione)
best_model <- glmnet(x, y, alpha = 1, lambda = best_lambda)

# Stampa i coefficienti del modello LASSO
coef(best_model)

## 11 x 1 sparse Matrix of class "dgCMatrix"
##                     s0
## (Intercept)  40.139921
## f1_z         -3.850296
## f2_z        -24.088259
## f3_z        -24.780290
## f4_z         19.280082
## f5_z        -13.008981
## f6_z          6.166835
## f7_z         14.584547
## f8_z        -10.168556
## f9_z          3.876559
## f10_z         5.714035

# Stampa il summary del modello lineare normalizzato precedente per confronto
summary(linModNormalized_1)

## 
## Call:
## lm(formula = won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + 
##     f7_z + f8_z + f9_z + f10_z, data = train_1)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -9.0156 -2.9560  0.1958  3.0363  9.3808 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  40.1581     0.1742 230.573  < 2e-16 ***
## f1_z        -18.9076     4.7109  -4.014 6.77e-05 ***
## f2_z        -48.9506     7.4952  -6.531 1.44e-10 ***
## f3_z        -27.7666     2.1909 -12.674  < 2e-16 ***
## f4_z         22.5459     1.9693  11.449  < 2e-16 ***
## f5_z        -13.0032     0.3099 -41.956  < 2e-16 ***
## f6_z          7.1655     0.7371   9.721  < 2e-16 ***
## f7_z         36.0699     6.4990   5.550 4.37e-08 ***
## f8_z         -9.8409     2.1973  -4.479 9.07e-06 ***
## f9_z          2.5699     1.2670   2.028    0.043 *  
## f10_z         7.4600     1.3101   5.694 1.98e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.206 on 573 degrees of freedom
## Multiple R-squared:  0.884,  Adjusted R-squared:  0.882 
## F-statistic: 436.7 on 10 and 573 DF,  p-value: < 2.2e-16

# Questo codice esegue stime usando un modello Lasso e un modello lineare (LM) su un set di dati di test e ne valuta le prestazioni utilizzando l'errore quadratico medio (RMSE) e le metriche    di classificazione binaria, ad esempio matrice di confusione, accuratezza, precisione, sensibilità, punteggio F e specificità.

# Predizioni utilizzando il modello Lasso
new <- data.matrix(test[, c("f1_z", "f2_z", "f3_z", "f4_z", "f5_z", "f6_z", "f7_z", "f8_z", "f9_z", "f10_z")])
prevLasso <- predict(best_model, s = best_lambda, newx = new)

# Predizioni utilizzando il modello LM (Linear Model)
new <- subset(test, select = c("f1_z", "f2_z", "f3_z", "f4_z", "f5_z", "f6_z", "f7_z", "f8_z", "f9_z", "f10_z"))
prevLM <- predict(linModNormalized_1, newdata = new)

# Calcolo dell'Errore Quadratico Medio (RMSE) per il modello Lasso
rmsLasso <- sqrt(mean((test$won - prevLasso)^2))

# Calcolo dell'Errore Quadratico Medio (RMSE) per il modello LM (Linear Model)
rmsLM <- sqrt(mean((test$won - prevLM)^2))

# Classificazione binaria utilizzando una soglia (0.5) per le predizioni del modello Lasso
prev <- ifelse(prevLasso > 0.5, "1", "0")
prev <- as.factor(prev)

# Matrice di Confusione per le predizioni del modello Lasso
confMatrix <- table(prev, test$won)

# Metriche di Performance
accuracy <- sum(confMatrix[1], confMatrix[4]) / sum(confMatrix[1:4])
precision <- confMatrix[4] / sum(confMatrix[4], confMatrix[2])
sensitivity <- confMatrix[4] / sum(confMatrix[4], confMatrix[3])
fscore <- (2 * (sensitivity * precision)) / (sensitivity + precision)
specificity <- confMatrix[1] / sum(confMatrix[1], confMatrix[2])

# Visualizza i Risultati
print(paste("Accuratezza:", round(accuracy, digits = 4)))

## [1] "Accuratezza: 0.6"

print(paste("Precisione:", round(precision, digits = 4)))

## [1] "Precisione: 0.6667"

print(paste("Sensibilità:", round(sensitivity, digits = 4)))

## [1] "Sensibilità: 0.6667"

print(paste("F-Score:", round(fscore, digits = 4)))

## [1] "F-Score: 0.6667"

print(paste("Specificità:", round(specificity, digits = 4)))

## [1] "Specificità: 0.5"

# Confronto tra RMS
print(paste("RMS Lasso:", round(rmsLasso, digits = 4)))

## [1] "RMS Lasso: 5.1466"

print(paste("RMS Modello Lineare:", round(rmsLM, digits = 4)))

## [1] "RMS Modello Lineare: 5.2047"

#TEST ANOVA

# Esecuzione di un test ANOVA per la variabile 'confID' come predittore
# Supponiamo che 'won' sia la variabile dipendente
resp_conf <- anova(lm(won ~ confID, data = train_1))

# Verifica se la variabile 'confID' ha un effetto significativo sulle vittorie
if (resp_conf["confID", "Pr(>F)"] < 0.05) {
  print("Si rifiuta l'ipotesi nulla perché la variabile Conference ha un effetto sulle vittorie. Posso inserire a modello la variabile.")
} else {
  print("Si accetta l'ipotesi nulla perché la variabile Conference non ha un effetto sulle vittorie. Non posso inserire a modello la variabile.")
}

## [1] "Si accetta l'ipotesi nulla perché la variabile Conference non ha un effetto sulle vittorie. Non posso inserire a modello la variabile."

# Esecuzione di un test ANOVA per la variabile 'divID' come predittore
# Supponiamo che 'won' sia la variabile dipendente
resp_div <- anova(lm(won ~ divID, data = train_1))

# Verifica se la variabile 'divID' ha un effetto significativo sulle vittorie
if (resp_div["divID", "Pr(>F)"] < 0.05) {
  print("Si rifiuta l'ipotesi nulla perché la variabile Division ha un effetto sulle vittorie. Posso inserire a modello la variabile.")
  div <- TRUE
} else {
  print("Si accetta l'ipotesi nulla perché la variabile Division non ha un effetto sulle vittorie. Non posso inserire a modello la variabile.")
  div <- FALSE
}

## [1] "Si accetta l'ipotesi nulla perché la variabile Division non ha un effetto sulle vittorie. Non posso inserire a modello la variabile."

#EFFETTI INTERAZIONE

# Verifica del valore di 'div' per costruire il modello appropriato
if (div) {
  linModNormalized1 <- lm(won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + f7_z + f8_z + f9_z + f10_z + divID + f1_z:f2_z + f4_z:divID, data = train_1)
} else {
  linModNormalized1 <- lm(won ~ f1_z + f2_z + f3_z + f4_z + f5_z + f6_z + f7_z + f8_z + f9_z + f10_z + f1_z:f2_z + f4_z:divID, data = train_1)
}

# Estrazione del riassunto del modello
riassunto <- summary(linModNormalized1)

# Estrazione dei nomi delle righe dal riassunto del modello
nomi_righe <- rownames(riassunto$coefficients)

# Nomi delle righe da escludere
nomi_da_escludere <- c("(Intercept)", "f1_z:f2_z", "f4_z:divIDCD", "f4_z:divIDMW", "f4_z:divIDNW", "f4_z:divIDPC", "f4_z:divIDSE", "f4_z:divIDSW", "divIDCD", "divIDMW", "divIDNW", "divIDPC", "divIDSE", "divIDSW")

# Creazione di un vettore con i nomi delle righe da utilizzare
nomi_righe_da_utilizzare <- setdiff(nomi_righe, nomi_da_escludere)

# Aggiunta della variabile 'divID' se presente nei nomi delle righe
if ("divIDCD" %in% nomi_righe) {
  nomi_righe_da_utilizzare <- append(nomi_righe_da_utilizzare, c("divID"))
}

# Indici di inizio e fine per le p-values della variabile 'divID'
div_indice_inizio <- which(rownames(riassunto$coefficients) == "f4_z:divIDCD")
div_indice_fine <- which(rownames(riassunto$coefficients) == "f4_z:divIDSW")

# Estrazione delle p-values della variabile 'divID'
div_p_values <- riassunto$coefficients[div_indice_inizio:div_indice_fine, "Pr(>|t|)"]

# Aggiunta di 'f4_z:divID' se la media delle p-values è inferiore a 0.05
if (mean(div_p_values) < 0.05) {
  variabili <- append(nomi_righe_da_utilizzare, c("f4_z:divID"))
}

# Aggiunta di 'f1_z:f2_z' se la sua p-value è inferiore a 0.05
if (riassunto$coefficients["f1_z:f2_z", "Pr(>|t|)"] < 0.05) {
  variabili <- append(nomi_righe_da_utilizzare, c("f1_z:f2_z"))
}

# Se nessuna delle condizioni precedenti è soddisfatta, utilizza solo le variabili esistenti
if (!(mean(div_p_values) < 0.05) && !(riassunto$coefficients["f1_z:f2_z", "Pr(>|t|)"] < 0.05)) {
  variabili <- nomi_righe_da_utilizzare
}

# Creazione della formula significativa per il modello
formula_significativa <- as.formula(paste("won ~", paste(variabili, collapse = " + ")))

# Costruzione del nuovo modello lineare con la formula significativa
linModNormalized1 <- lm(formula_significativa, data = train_1)

# Visualizzazione del riassunto del nuovo modello
summary(linModNormalized1)

## 
## Call:
## lm(formula = formula_significativa, data = train_1)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -9.6428 -2.8848  0.1359  2.9845  9.4655 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  39.9445     0.1972 202.572  < 2e-16 ***
## f1_z        -21.3191     4.8113  -4.431 1.12e-05 ***
## f2_z        -51.5660     7.5554  -6.825 2.25e-11 ***
## f3_z        -26.5996     2.2420 -11.864  < 2e-16 ***
## f4_z         21.6454     2.0014  10.815  < 2e-16 ***
## f5_z        -13.0055     0.3088 -42.117  < 2e-16 ***
## f6_z          6.6184     0.7726   8.566  < 2e-16 ***
## f7_z         38.5394     6.5652   5.870 7.38e-09 ***
## f8_z         -9.3337     2.2005  -4.242 2.59e-05 ***
## f9_z          2.5617     1.2624   2.029   0.0429 *  
## f10_z         6.4646     1.3764   4.697 3.31e-06 ***
## f1_z:f2_z    -0.4396     0.1927  -2.281   0.0229 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.191 on 572 degrees of freedom
## Multiple R-squared:  0.8851, Adjusted R-squared:  0.8829 
## F-statistic: 400.4 on 11 and 572 DF,  p-value: < 2.2e-16

Può aiutare l’uso della POISSON?

# Questo codice addestra un modello di regressione generalizzata di Poisson (linModNormalized2_pois) utilizzando la stessa specifica del modello lineare (linModNormalized2). Successivamente,    vengono ottenuti e combinati i coefficienti di entrambi i modelli in un unico dataframe per una facile comparazione.

# Creazioni di un modello di regressione generalizzata di Poisson
linModNormalized1_pois = glm(formula_significativa, family=poisson(link=log), data = train_1)
summary(linModNormalized1_pois)

## 
## Call:
## glm(formula = formula_significativa, family = poisson(link = log), 
##     data = train_1)
## 
## Coefficients:
##               Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  3.6487423  0.0076641 476.082  < 2e-16 ***
## f1_z        -0.0759700  0.1871492  -0.406 0.684792    
## f2_z        -0.6559339  0.2949496  -2.224 0.026156 *  
## f3_z        -0.8218722  0.0881466  -9.324  < 2e-16 ***
## f4_z         0.6768258  0.0787680   8.593  < 2e-16 ***
## f5_z        -0.3390626  0.0120410 -28.159  < 2e-16 ***
## f6_z         0.1840628  0.0281955   6.528 6.66e-11 ***
## f7_z         0.3814903  0.2581319   1.478 0.139437    
## f8_z        -0.2773063  0.0821204  -3.377 0.000733 ***
## f9_z         0.0838793  0.0479279   1.750 0.080098 .  
## f10_z        0.2164135  0.0520268   4.160 3.19e-05 ***
## f1_z:f2_z   -0.0005862  0.0072869  -0.080 0.935880    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 2294.39  on 583  degrees of freedom
## Residual deviance:  326.91  on 572  degrees of freedom
## AIC: 3552
## 
## Number of Fisher Scoring iterations: 4

# Ottenimento dei coefficienti per ambo i modelli
normal = coefficients(linModNormalized1)
poisson = exp(coefficients(linModNormalized1_pois))

# Combinazione dei coefficienti in un unico dataframe
coefficients_table <- cbind(normal, poisson)

coefficients_table

##                  normal    poisson
## (Intercept)  39.9445260 38.4263085
## f1_z        -21.3191022  0.9268440
## f2_z        -51.5659536  0.5189572
## f3_z        -26.5996116  0.4396079
## f4_z         21.6454468  1.9676222
## f5_z        -13.0055452  0.7124379
## f6_z          6.6184349  1.2020914
## f7_z         38.5393686  1.4644655
## f8_z         -9.3336530  0.7578224
## f9_z          2.5617476  1.0874976
## f10_z         6.4645993  1.2416157
## f1_z:f2_z    -0.4395592  0.9994139

#Considerazioni Generali:
# È importante notare che gli effetti delle variabili possono essere interpretati in modo diverso a seconda della distribuzione scelta per il modello. La scelta tra distribuzione normale e di   Poisson dipende dalla natura della tua variabile dipendente e dai tuoi obiettivi di modellazione. Nel contesto di modelli di regressione, è sempre buona pratica verificare l'adeguatezza del   modello esaminando i residui, eseguendo test diagnostici e valutando la bontà di adattamento. L'interpretazione dei coefficienti dovrebbe essere fatta considerando la scala appropriata per    la distribuzione utilizzata (lineare per la normale, logaritmica per la Poisson). Se stai cercando di prevedere il numero di vittorie, la distribuzione di Poisson potrebbe essere più
# appropriata per variabili conteggio come questa. Tuttavia, è sempre necessario verificare l'adeguatezza del modello ai dati specifici.